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ABSTRACT 



For a universe with massive neutrinos, cold dark matter, and baryons, we solve the linear 
perturbation equations analytically in the small-scale limit and find agreement with numerical 
codes at the 1 — 2% level. The inclusion of baryons, a cosmological constant, or spatial curvature 
reduces the small-scale power and tightens limits on the neutrino density from observations of 
high redshift objects. Using the asymptotic solution, we investigate neutrino infall into potential 
wells and show that it can be described on all scales by a growth function that depends on 
time, wavenumber, and cosmological parameters. The growth function may be used to scale the 
present-day transfer functions back in redshift. This allows us to construct the time-dependent 
transfer function for each species from a single master function that is independent of time, 
cosmological constant, and curvature. 



Subject headings: cosmology: theory - dark matter - large-scale structure of the universe 



Introduction 



The mixed-dark-matter (MDM) scenario for structure formation involves a hot component of massive 
neutrinos as well as the usual cold and baryonic dark matter components. In this case, even calculations in 
linear perturbation theory are non-trivial due to the time-dependent energy-momentum relation and non- 
vanishing angular moments of the neutrino distribution. Perturbations no longer grow uniformly with time 
independent of scale. Specifically, the growth of fluctuations is suppressed below the time-dependent free- 
streaming scale of the neutrinos due to collisionlcss damping. Numerical calculations with state-of-the-art 



Boltzmann codes (e.g. Ma & Bertschingei 1995; Dodelson, Gates, & Stebbins 1996; Seljak & Zaldarriaga 



1996) still require a fair amount of time to solve the evolution equations on small spatial scales. Moreover, 
the additional parameter represented by the neutrino mass m„ makes an exhaustive search of parameter 
space more difficult and has led most workers to date to fix parameters such as the baryon density (see e.g. 
Ma 1996). For these reasons, we consider here an analytic treatment of small-scale perturbation theory in 
MDM cosmologies. 

The inclusion of baryonic dark matter further complicates the dynamics. Recent measurements of 



high-redshift deuterium abundances (Tytler et al. 1996; but see Rugers fc Hogan 1996) and new theoretical 



interpretations of the Lyman-a forest ( Weinberg et al. 1997, and references therein) suggest a value of the 
baryon density fib greater than the fiducial nucleosynthesis value of 0.0125/i~ 2 ( Walker et al. 1991). Baryons 
suppress fluctuations on small scales because, prior to recombination, photon pressure from the cosmic 
microwave background supports them against collapse. Hu & Sugiyama (1996, hereafter |HS96j ) developed 
a formalism to account for this effect and solve the evolution equations exactly on small scales. The key 
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aspect to the treatment is the ability to ignore completely the role of baryons as a gravitational source for 
enhancing CDM fluctuations. 

Massive neutrinos act in a manner similar to the pre-recombination baryons. On scales smaller than 
their free-streaming length, the neutrinos are smoothly distributed and hence do not contribute to the growth 



of perturbations. Here we generalize the techniques of HS96 to include the hot component, thereby allowing 
us to solve analytically for the transfer function on the smallest scales. We then consider how to describe the 
end of free streaming and the resulting infall of neutrinos into the existing potential wells. This allows us 
to collapse all of the late-time neutrino effects and base the transfer function on a single time-independent 
function of scale. 

In an MDM cosmology with realistic baryon content, the amplitude of small-scale fluctuations is impor- 
tant due to growing evidence of early structure formation from high-redshift observations. The model has 
difficulty in explaining observations of galaxies at redshift z ~ 3 (Steidel et al, 



as well as damped Lyman-a systems at at comparable redshift (Mo & Miralda-Escude 1994; Kauffmann 



1996; 



Mo fc Fukugita 1996) 



fc Chariot 1994; Klypin et al. 1995; Ma et ak| 1997). Baryons only exacerbate this problem and tighten 
the upper limit on f)„. Indeed, they yield a stronger effect for MDM as compared to CDM cosmologies 
because Oj, / (Qo — ®>v) rather than f^/f^o enters into the fluctuation amplitude. Similarly, the growth rate of 
fluctuations is determined by Q^/D,q such that a given £!„ causes more suppression in a low-density universe. 
Our results here should therefore aid in the investigation of the parameter space left available to MDM 
cosmologies. 

The outline of this paper is as follows. After establishing the notation in §||, we present in 

§1 the 

small-scale solutions of the perturbation equations derived in the Appendix. We use these solutions in §^ 
to study the behavior of neutrino infall and to find analytic approximations thereof. From these results, we 
construct in §[s] the transfer functions in time and wavenumber for the cold dark matter and total density 
perturbations and find agreement at the percent level with analytic results in the small-scale limit. In §|^, 
we show how these results may be scaled to cosmologies with a cosmological constant or spatial curvature. 



2. Notation 

We begin by establishing the notation used throughout. The density of the ith particle species (i = c, 
cold dark matter; b, baryonic dark matter; massive neutrinos) today in units of the critical density is 
denoted f2j, whereas the fraction of the total matter density today f2o = Y)j &j is denoted /$ = fii/fio- 
As short-hand, we employ for example i — cb to denote f c b = f c + fb- Note that f c + fb + fv = 1- 
Density perturbations are expressed as Spi/pi = Si, where the hybrid combinations are density weighted 
(e.g. 5 c b = fc.S c + fbSb)- The CMB temperature is given by Tqmb = 2.7O2.7K; the best determination to 



date is 2.728 ± 0.004 K (Fixsen et al. 1996; 95% confidence interval) at which it is fixed for most of our 
expressions. Finally, as usual the Hubble constant is written as Ho = lOOh kms -1 Mpc -1 . 

Time is parameterized as y = (1 + z cq )/(l + z), where 

z cq = 2.50 x io A n a h 2 e^ (i) 



is the redshift of matter-radiation equality. The second important epoch is when the baryons are released 
from the Compton drag of the photons near recombination, i.e. yd — y(zd) where (HS96; Eiscnstcin & Hu 
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1997a) 

/o 1^0.251 

- = 12 ^ + o: ( ^r |1+tl( ^ )feil - (2) 

6i = 0.313(r!o^ 2 )~ - 419 [l + 0.607(r! / l 2 ) - 674 ], 
b 2 = 0.238(fi /i 2 ) - 223 . 

After this epoch, baryons fall into the potential wells provided by the cold dark matter and participate in 
gravitational collapse. 

We often label the comoving wavenumber k relative to the scale that crosses the horizon at matter- 
radiation equality, thus defining the quantity 

e 2 .,^ 2 )- 1 = -^-(fio#o 2 r 1/2 (i + ^r 1 ' 2 ■ (3) 



Mpc- 1 z,v ' 19.0 

The small-scale limit is defined as q 3> 1. In the next section (^), we place an additional restriction that 
the momentum of the neutrinos keep them out of the perturbations formed by the heavier species. Such 
scales are below the free-streaming scale, which itself shrinks with time [c.f. eq. We show in §Q how to 

account for neutrino infall. 

We often encounter functions of wavenumber or time that depend additionally on cosmological param- 
eters; we denote these as e.g. F(y, q; f u . . .). Where the cosmological parameter dependence is not being 
emphasized, we often drop the parameters after the semicolon, e.g. F(y, q) — F(y, q; j v . . .). 



Small Scale Solution 



Below the free-streaming scale of the neutrinos (see §[|) and sound horizon of the baryons at recom- 
bination, the equations of motion for matter density fluctuations may be solved analytically in a matter 

The key approximation is that on sufficiently small 



radiation universe using the techniques of HS96 



scales the neutrinos move too quickly to trace the perturbations in the CDM and baryons. In this case, the 
neutrinos contribute no gravitational sources to the evolution equations of the other species, thereby slowing 
the growth of fluctuations. The baryons have a similar behavior prior to the drag epoch; in the Appendix 
we describe how to include both effects. The result is that density perturbations grow as 



3cb(y, q; U, fb, yd) = D cb (y; f v )5 d {q\ /„, f b , y d ) 



(4) 



Equation (||) states that the density perturbation today is the product of a growth function D c b that depends 
on the neutrino fraction and the amplitude 5^ of fluctuations entering the growing mode at the Compton 
drag epoch yd- 



The quantity 



Pi(fi) 



VI + 24/, 



> 



(5) 



determines the reduction from a linear growth rate S oc y, where fi is the fractional density in gravitationally 
clustering components: i — c and i = cb before and after the drag epoch respectively. Hence the growth 



factor is given by |Bond et alj (1980) 



Dcb{y-Jy) = y 



l-Pc 



(6) 
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where we take fio = 1; we generalize to fio ^ 1 in §g. 

The amplitude of the fluctuation entering the growing mode at the drag epoch is 

$d(q; f v , fb, Vd) = 9.50M d (g)$(0, q) , 



(7) 



where $(0, q) is the initial amplitude of the potential perturbation. The quantity Md expresses the matching 
condition between the growing mode y 1 ~ Pc before the drag epoch and y 1 ~P^ after the drag epoch, as well 
as the matching required to describe the onset of matter domination. We find 



M d (q; fv, fb,Vd) 



f c 5 - 2(p c +p cb ) 
f cb 5 - Ap cb 



(i + y d y° b - p ° 



1 



(3-4p c )(7-4p c6 ) 



(i + ydY 



Ai(q). 



where0 



with 



M{q\Ub,fu) 



l- 0.553/„ b + 0.126/3 
— ln(1.84/3 c(? ) 

1-0.193/*^ + 0.169/„ 



(8) 



(9) 



& = [1 - Q.949/^]- 1 . (10) 

Equation (||) results from a series expansion in (1 + L/d) 1 of the analytic solution and accounts for small 
deviations from the power-law growing-mode behavior due to radiation. The expansion is only accurate for 
fb + fv < 0.6 and y d ;> 1, but the form in the Appendix is general. Notice that Md — > A\ as fb — > and 
that the term in brackets introduces only a small correction since y d p> 1 in cases of interest. 

In Fig. |], we display an example of the time evolution of a mode given by the analytic solution (including 
here the decaying mode given in the Appendix, important for y < y d ) compared with numerical solutions. 
The 5% offset at early times is due to changes in the expansion rate as the neutrinos become non-relativistic. 
The oscillatory errors arise from neglect of baryon acoustic oscillations; these Silk damp away well before the 
drag epoch for these scales. The main effect of the neutrinos is to slow the growth of the CDM and baryons 
after equality {y > 1) because they represent a smooth gravitationally-stable component on these scales. 
Similarly, since the baryons have no fluctuations on these scales until they fall into CDM potential wells 
after the drag epoch yd, they reduce the growth rate between equality and the drag epoch. Furthermore, 
they reduce the net fluctuation S c b as f c / fcb leading to the offset between the curves in Fig. [|. 



4. Neutrino Infall 



Eventually, the neutrinos fall into the CDM potential wells, breaking the approximation of the last 
section. The neutrino thermal velocity decays with the expansion of the universe as v v cx (am^)^ 1 : infall 



occurs when their velocity slows sufficiently as to allow clustering by the Jeans criteria (Bond & Szalay 1983; 



Ma 1996): 



fcf s = 4 7 rGpa>*oc(l + z)- 1 / 2 (^) 3 / 2 (/,/Ag 
qis « y 1/2 (U/N„). 



(11) 



2 We assume here the number of massive neutrinos 7V„ = 1; see equation (A17) for the general case 
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Fig. 1. — Growth suppression from the neutrinos and baryons. The addition of neutrinos slows the growth of 
CDM + baryon fluctuations after matter-radiation equality y > 1 (upper panel, curve b compared with curve 
b). Increasing the baryon fraction further suppresses fluctuations by a time- independent factor at y > yd 
due to growth suppression in the CDM fluctuations for 1 < y < yd and the lack of baryon fluctuations in the 
weighted density perturbations (upper panel, curve c). The analytic expressions (dashed lines) agree with 
numerical results (solid) for Oo = 1, h = 0.5 and q = 160 at the 1% level except for 5% discrepancies at 

y<i. 



Recall that k and q are related by equation (^). Here N v is the number of massive neutrino species, assumed 
to be degenerate in mass. For simplicity, we restrict our examples to N v = 1 throughout, but we have 
verified that the infall description is valid for N v ^ 1 at the 1 — 2% level in the range < z < 25.^ 

On scales q -C q{ s , the neutrinos will follow the cold dark matter. By acquiring density perturbations, 
they enhance the CDM+baryon potential wells and drive the growth rate back up to y. The CDM+ baryon 
density fluctuations thus acquire a scale dependence to their growth rate 



D C b{y,q\f v ) 



1 



y 



l 



Pcb/0.7 



yts(q;fu) 

where the coefficient 0.7 represents a fit to the numerical evolution 



V 



1— Pcb 



(12) 



The transition epoch j/f s incorporates two effects. The first is that from equation (|Tl]) the characteristic 
epoch for infall for a given wavenumber must scales as y oc {qN v / f„) 2 . The second is that the growing modes 



3 Effects at redshifts approaching the epoch at which the neutrinos become non-relativistic are not accounted for by this 
approximation. For N v ^ 1, the growth function of equation ( |l2| ) allows scaling at low redshifts but this does not imply 5^ or 
Tmaster is independent of 7V„ around the maximal infall scale. We treat these effects in Eisenstein & Hu (1997b). 
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Fig. 2. — Analytic vs. numerical descriptions of neutrino infall for f u — 0.3, fb — 0.01, Slo = 1, h = 0.5. Upper 
panel: Growth functions compared to the fully free-streaming time-dependence of y 1- ^; numerical CDM 
+ baryon fluctuations compared with growth function D c b (short- dashed lines) numerical CDM + baryon + 
neutrino density- weighted fluctuations compared with D c b v (long-dashed lines). No infall is represented by 
a horizontal line of amplitude 1 (cb) and f c b — 0.69 (cbv) here. Lower panel: relative error. 

of the free-streaming epoch oc y 1 ^P<='> and the infall epoch oc y must be matched across the transition. This 
matching condition may only depend on /„. We find that the total effect is well approximated by 

y {s (q; fu) = 17.2/„(1 + 0A88f- 7 ^)(qN„/ j v f , (13) 

In Fig. |^ we show the enhancement of density fluctuations due to infall and compare D c b from equation ( |l^ ) 
to numerical solutions for various q (short-dashed lines). 

Furthermore, the density-weighted fluctuation 

Scbu = fcbdcb + fy5 y (14) 

follows from the infall solution by noting that it converges above the free-streaming scale to S c b and below 
to fcbScb- Thus 

<W(y, q; h, yd) = D cbv (y, q; fv)8d(q\ fv, h, yd) , (15) 

where 

Pcb/0.7 

y 1 -**. (16) 
In Fig. 0, we show that this form produces a good fit to the numerical results (lower set, long-dashed lines). 



D cbv (y,q;f v ) 



f0.7/pcb 
Jcb 



0.7 



v l + Vf s (q;f v ) J 
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Fig. 3. — Demonstration of the existence of a time-independent master function. The transfer functions 
for the CDM + baryon (T c b) and CDM + baryon + neutrino (T c b v ) fluctuations at 3 different rcdshifts are 
divided by the growth factors D c b and D c b u respectively to obtain estimates of T mas t C r. That the 6 curves 
superimposed in the upper panel agree at the 1 — 2% level (relative to T c b/D c b at z = 0, as shown in the lower 
panel) establishes the existence of the master function and verifies the accuracy of the growth functions. The 
analytic prediction for T mastcr (long-dashed line) converges to within 1% of the numerical results at q ^> 1. 
The model here is Q = 1, h = 0.5, /„ = 0.4, f b = 0.2. 



Finally, the horizon at the epoch when the neutrinos become nonrelativistic sets the maximal free- 
streaming scale. Beyond this scale, the neutrinos are always in the infall regime and the evolution of density 
fluctuations becomes independent of the neutrino fraction 

limS cb (y,q) = Y ^ym q ). (17) 

The appearance of y/(l + yi s ) in equation ( |l2] ) assures the proper time dependence for the evolution in the 
large-scale limit. The growth functions are thus not subject to a small-scale approximation and remain valid 
for all q. 



5. Transfer Functions 

We are now in a position to evaluate the transfer function, defined as 
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for the ith component of the matter. For the CDM + baryon (i — cb) and the CDM + baryon + neutrino 
(i = cbv) systems, we obtain 

Tcb(y,q;fv,fb,ijd) = y~ 1 D cb (y,q;f v )T maBta: (q- > f t/ ,f b ,ya), 
Tcbv(y,q;fv,fb,yd) = y~ 1 D cbl/ (y,q;f l/ )T mastcI (q;f l/ ,f b ,y d ), (19) 

where T mastor has the small-scale limit of 



Urn T mastor (q; f v ,f b , y d ) = - - - 2 , (20) 



M d (g) 

q^oo 14.4^ 

which follows from equations (Q) and (|l7j). It is important to note that relation ([H]) holds independently of the 
small-scale approximation: namely that once the growth factors y~ x Di are removed, the transfer functions 
depend only on a single time-independent function of q related to perturbations in the CDM component 
at the drag epoch. This simplification holds only for y 3> yd] near the drag epoch, the contribution of the 



decaying mode cannot be neglected. In Eiscnstein & Hu (1997b), we exploit the existence of this master 
function to obtain the full time- and g-dependent transfer functions based on fitting formulae for T master . 

We show a comparison between analytic and numerical results in Fig. ||. The numerical T cb , T cbv 
(lower curves, solid and dashed lines respectively) at 3 different redshifts are plotted here. The upper curve 
represents the master function T mas t cr obtained through inverting the growth factor. Note that 6 different 
estimates of T mastC r obtained from T cb and T cbll are superimposed and agree at the ~ 1 — 2% level. Also 
shown is the small-scale prediction of equation (120) which converges rapidly for q J> 1. 



Finally, note that the neutrino transfer function is implicitly defined as 

Tu = f- 1 {T cbv -f cb T cb ), (21) 

and its construction from the growth functions and T master yields density-weighted errors on the same order 
as the other transfer functions, i.e. (1 — 2%) x 8 cbv / f v 5 v . 



6. Low Density Models 



The formulae presented thus far are valid for 0,^ — 1 cosmologies. We have however shown that the 
transfer functions today can be expressed as the products of a growth function and a master function 
related to fluctuations at the drag epoch. At this earlier time, the cosmological constant and curvature have 
negligible effects on the dynamics. This implies that a simple modification of the growth function at late 
times to account for SIq 7^ 1 effects will suffice for a complete description. 

Let us recall that on the largest scales, where neutrino free streaming and radiation pressure gradients 



are negligible, fluctuations grow as ( Heath 1977) 



5 f v dx 
D(y,n ,tt A ) = -g(y) J 



(22) 



where 



9 2 {y) = 2T 3 + 2T 2 % Hi - «o - ^a)/U + Vo 3 ^a/^o , (23) 
with yo — (1 + z C q). Analytic forms for the J7a = and f2o + = I cases are given in proth & Peebles] 



(1975) and Bildhauer et al. (1992) respectively. The normalization of the growth rate has been chosen so 
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that D — y at early times. Moreover, equation (|22j) states that after matter ceases to dominate the expansion 
rate, fluctuation growth halts. 

By matching these asymptotic limits, we can approximate the growth function in the presence of neu- 
trinos by the replacement of y with D, i.e. 



D C b{y,q; fv,tt ,n A ) 



l + 



D{y) 



i + «&(«) 



0.7 



Pcb/0.7 



f 0.7/pcb 
Jcb 



D(y) 



0.7" 



Pat/0.7 

A + yis(q)J 

where of course D(y) = D(y; Qo,fl\) implicitly. Likewise, the transfer functions become 

T c b(y,q) = D~ 1 (y)D cb (y,q)T masta (q) , 
T c bv{y,q) = D- 1 (y)D cb „(y,q)T mastel (q) . 



(24) 
(25) 



(26) 



We show an example in Fig. [|. The model has been chosen to have Qq = 0.25, Da = 0.75 with the same 
Qoh 2 i fb and /„ as Fig. [3] and hence the same Tmastcr- The invariance of the master function is demonstrated 
by overplotting estimates from T cb and T cbv in this model and from T cb of Fig. || The lower panel shows the 
fractional difference of the former with respect to the latter. 

In principle, the fact that infall freezes out at a later redshift in a A vs. open cosmology of the same 
fio allows one to distinguish between these two alternatives by the shape of the transfer function at a single 
redshift alone. Massive neutrinos therefore break the shape degeneracy of the transfer function in low-density 
universes. However, since this effect is only significant if Qq <C 1, the more important fact is that growth in 
the matter-dominated epoch depends on (1 — VL V / S7o) and hence even a relatively small density in neutrinos 
can make a dramatic effect on the growth rates in a low-density universe. Fig. || demonstrates that il„ = 0.1 
in an Oq = 0.25 universe gives the same magnitude effects as f2„ = 0.4 in an f2g = 1 universe. Upper limits 



on VL V from small-scale fluctuations thus tighten for low density universes (Primack et al. 1995) 



7. Discussion 

We have presented small-scale solutions to the linear evolution of perturbations in a cold + hot + 
baryonic dark matter cosmology. They converge to within 1% of the numerical solutions for 5 » 1 and 
q 3> <?fs and are expressed in terms of hypergeometric functions in the Appendix. We have also given 
simplified forms in the main text, employing only algebraic functions, which are valid for + f2„ < O.6f2o 
and floh 2 > 0.1. Baryons play a significant role in MDM cosmologies since dark-matter perturbations are 
density weighted and depend on f2b/(f2o — f2„) not f2f,/f2o- Likewise neutrinos in a low-density universe yield 
enhanced effects since growth rates depend on f2„/f2o. 

By comparing analytic and numerical results at q <J q{ s , we have isolated the effects of neutrino infall 
and described them to an accuracy of 1 — 2% by time- and wavenumber-dependent growth factors for the 
CDM and total matter. The freeze-out of infall at late times when A or curvature come to dominate can 
similarly be taken into account. We have shown that these growth factors are valid beyond the small-scale 
approximation. 
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Fig. 4. — Demonstration of the invariance of the master function (T maste r) under changes in the cosmological 
constant A. The model has the same Oo/i 2 , /& and f v as Fig. [| at z = but with f2o = 0.25. Division of T c t, 
and T c bu by the growth functions returns the master function T mas ter over which the estimate of Fig. [5] and 
the small-scale analytic solution are plotted (top panel). The errors are plotted relative to T mas t e r from the 
T c fc, z — estimate of Fig. || and show agreement at the 1 — 2% level (bottom panel). 



The full time and wavenumber dependent transfer functions for the CDM and total matter can thus be 
described as a product of these growth factors and a single function of wavenumber related to the amplitude 
of fluctuations at the Compton drag epoch that in turn depends on floh 2 , ^W^Oj ^b/^Oi an d the number 
of degenerate neutrino species. We leave a quantitative description of this master function and implications 
for constraints on f2„ to a companion piece flEiscnstcin fe Hu| 1997b). 

Acknowledgments: W.H. and D.J.E. are supported by NSF PHY-9513835. W.H. was additionally supported 
by the W.M. Keck Foundation. Numerical results were extracted from the CMBfast package of Seljak fe 



Zaldarriaga (1996) v. 2.3. 



Derivation of the Small Scale Solution 



Following the analytic approach of HS9(: , one can solve the equation of motion for the CDM on small 
scales where the gravitational effects of the MDM can be neglected. The idea is to separate the evolution 
before and after the drag epoch. Before the drag epoch, the gravitational effects of the baryons can be 
ignored below the sound horizon as they are pressure supported by the photons. The equation of motion 
then becomes 

6c + -Sc = ^fctooH%8 c , (Al) 
a la 
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where ovcrdots indicate derivatives with respect to conformal time rj. Unfortunately, the complicated equa- 
tion of state of massive neutrinos prevents the time evolution of the scale factor a(r)) from being written down 
in closed form. However, we know that the neutrinos behave as radiation at early times when T v 3> m v 
and as matter at late times when T u -C m v . These limits are identical to those found if one considered 
the neutrinos to be massless and added their mass to that of the non-relativistic matter. It is therefore a 
reasonable first approximation to leave the background evolution unmodified by the replacement of a portion 
of the non-relativistic matter with massive neutrinos, i.e. 

77 = 2(fi # 2 r 1/2 (i + W 1/2 (Vi + w - 1) > (A2) 

where y = (1 + z cq )/(l + z). Here, z cq [eq. (Q)] assumes three massless neutrino species with the usual 
thermal history. We neglect cosmological constant and curvature effects here (see §^). This form otherwise 
errs only between z eq and the epoch at which the neutrinos become non-relativistic. Our approach will be to 
use this approximation to solve the small-scale limit analytically and then correct for a <J 5% modification 
due to changes in the expansion rate. 



With this approximation, equation (Al) can be rewritten in terms of y = (1 + z cq )/(l + z) as 

d\ (2 + 3y) d 3 
dy 2 2y(l + y)dy 2y(l+y) 

As shown in ES96| , the general solution to this equation is given through Gauss's hypergeometric function 
F (also written 2-F1) 

Ui = {l + yy^Ffaoi + \,2a t + i; -L-) , (A4) 

2 2 1 + y 



where i = 1, 2 and 

1 T VI + 24/ c 



4 (A5) 

with — and + for i = 1 and 2, respectively. It is useful to note that liniy^oo Ui — y~ ai and thus the two 
solutions represent the growing and decaying mode for i = 1,2 respectively. Clearly, ol\ + a-z = 1/2. 

We obtain the amplitude of the growing and decaying mode by matching onto the y -C 1 solution of 
HS96| eq. (B12)| 

5 C = 9.501n(9.24,7y)$(0, g ) , y < 1 , (A6) 

to find 

5 c (y,q) = 9.50[A 1 (q)U 1 (y) + A 2 (q)U 2 (y)M0,q) y<y d , (A7) 
where the matching coefficients are 

A ^ = r(o T{ 7uo^ + 1 !S \ M$Mq) + 2V,(1) - 2^(1 - 2 ttl ) - 2 In 2] , (A8) 
I (2ai + l/2)27rcot(27rai) 



with ip( x ) = r'(x)/r(x). The expression for Ai follows from equation (AS) with the replacement 1 — > 2 in 
the subscripts. 



4 The numerical factors here reflect the kick a perturbation gets at horizon crossing deep in the radiation-dominated epoch. 
They are cali brated to agree at the 1% l evel with C MBfast v. 2.3 (high precision version) and represent a 1-2% shift vers us the 



calibration of Hu et al 
at the 1% level. 



(1995) based on Bugiyama (1995). Our calibration also matches the code of M. White (Hu et al. 1995) 
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Equation (A7) takes the evolution up to the drag epoch. After the drag epoch, the baryons are released 
from the photons and behave as cold dark matter. The equation of motion for the combined CDM+baryon 
system is the same as equation dA3) but with the replacements 



5 C — * -f-^c + -^-S b = 5 cb . 

Jcb Jcb 



fc 



cb ■ 



and therefore has the same solutions as given by equation (A4) with arguments 



1 ± VI + 24/, 



cb 



(A9) 



(A10) 



The full solution for 5 c b is now found by matching the solution of equation (A7) onto the new growing 
and decaying modes. Since the decaying mode becomes unimportant for y 3> y d , 



where 



5 cb (y,q) = 9.50M d (q)U 1 (y)$(0 > q), 



d Q fcb\ IW2 - U[U 2 



(All) 
(A12) 



gives the matching condition. 



Furthermore, the baryons fall into CDM potential wells at the drag epoch and subsequently follow them 
so that 

Scb = S c = S b , z/ > y d . (A13) 



For high fio/1 2 ;> 0.1 (or h ;> 0.3), y d ^S> 1 so that the decaying-mode contribution to equation (All) 
coming through U2 is negligible. Furthermore the growing mode can be expanded in powers of (1 + y d ) _1 as 



M d (q) 



f c 2(5i + <*!)-! 

jcb ^OL\ - 1 



(i + w) ai - 



l + ii + Vd)- 1 ^^ 1 1 



1 



(4ai + l)(45x -3) 



A x {q). 



(A14) 



The approximation holds for f c J> 0.2. Note that the apparent divergence at f c = 0.125 or a,\ = —1/4 
disappears once the decaying mode is properly included. 



The quantity A\ can either be evaluated from equation (A8) or approximated as 

Ax{q) = [1 - 0.553/U + 0.126/, 3 b ] ln(1.84/3 c g) , 



(A15) 



with 



I3 C = exp[-2V>(2a 2 ) - 2^(3)] « [1 - 0.949/^r 1 , (A16) 

which is valid at the 1% level for f vb < 2/3. Finally, we correct for the change in the background expansion 
rate by comparing equation (All) to an numerical integration of equation (Al) and find the small correction 

AM 



Ax(q) 



1 -0.193V7PV^ + 0.169/^°- 2 ' 



(A17) 



for f v < 0.6. With the identities a% = p c — 1 and a\ = p c b — 1, the derivation of the small-scale evolution 
[eq. (||)] is now complete. 
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